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Abstract: In this work we study the process of mutant invasion on an example of a 
consumer- resource system with annual character of the behavior. Namely, individuals are 
active during seasons of fixed length separated by winter periods. All individuals die at 
the end of the season and the size of the next generation is determined by the number 
of offspring produced during the past season. The rate at which the consumers produce 
immature offspring depends on their internal energy which can be increased by feeding. The 
reproduction of the resource simply occurs at a constant rate. 

At the beginning, we consider a population of consumers maximizing their common 
fitness, all consumers being individuals having the same goal function and acting for the 
common good. We suppose that a small fraction of the consumer population may appear 
at the beginning of one season and start to behave as mutants in the main population. We 
study how such invasion occurs. 
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Se rallier ou filouter ? 
Une analyse evolutionnaire dans un systeme 
consommateurs-ressources. 

Resume : Dans ce travail nous etudions les processus evolutifs a l'oeuvre dans un sys- 
teme consommateurs-ressources "saisonnier" dans lequel les consommateurs ont un compro- 
mis dynamique a faire entre le temps alloue a la recherche de nourriture et celui alloue a 
leur reproduction. Les individus sont actifs pendant des saisons de longueur fixe, separees 
par des periodes d'hiver ou seuls les immatures produits durant la saison survivent (oeufs, 
graines,...). La taille de la generation suivante de matures est alors determinee par ce nombre 
de survivants. 

Dans ce rapport, nous considerons d'abord une population de consommateurs qui maxi- 
misent leur "fitness" commune (l'efficacite reproductive de la population entiere), tous les 
consommateurs ayant ce meme objectif; ils agissent en quelque sorte pour le bien commun. 
Nous supposons par la suite qu'une petite fraction de la population des consommateurs peut 
apparaitre au debut d'une saison et choisir une strategic differente de celle de la population 
principale; nous appellerons ces individus deviants "mutants". Nous etudions en details la 
strategie mise en oeuvre par ces mutants et analysons pour finir leur capacite ou non a 
supplanter la population residente. 

Mots-cles : Commande optimale, Jeux differentiels, Modele consommateur-resource, 
Invasion 
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1 Introduction 

Biodiversity found on Earth consists of millions of biological species, thousands of different 
ecosystems. Among this variety, one can easily identify many examples of resource-consumer 
systems like prey-predator/parasitoid-host systems known in biology |T| or birth-death sys- 
tems known in cell biochemistry [2J. Usually, individuals involved in such kind of systems 
(bacteria, plants, insects, animals) have conflicting interests and models describing such in- 
teractions are based on principles of game theory [21 IH [SJ |5] . Hence the investigation of these 
models is of interest both for game theoreticians and for biologists working in behavioral 
and evolutionary ecology. 

One of the first questions that usually appears when first consulting evolutionary the- 
ory books is the following: could we say that individuals behave rationally or optimally 
throughout their life? The answer is most probably "yes" if we consider the evolution as a 
slow process tending to some equilibrium. Following Darwin theory and its main statement 
about the survival of the fittest, we can assume that evolution of populations leads to a 
situation where individuals maximize their fitness or try to protect themselves from invasion 
by others [7: 8J. Such population can be referred to as residents who use an optimal maxi- 
mizing strategy or an uninvadable strategy respectively. The first type of strategy could be 
dynamically stable and lead to an asymptotically stable equilibrium, but it could also not 
be. Particularly, and this is a well-known fact in economics, a free-rider may overcompete 
competitors cooperating with him by "cheating" and using a "greedy" strategy. In the se- 
quel, populations which behave differently from the residents will be termed mutants. On 
the other hand, if residents use an evolutionary stable strategy, this will not allow them to 
get the maximum possible value of the fitness but will help them avoid mutant invasion. 
This seems reasonable from a biological point of view but such strategy could be dynamically 
unreachable or could not lead to a stable equilibrium in a long-term perspective [3] . 

In this work we study the process of mutant invasion on an example of a consumer- 
resource system with annual character of the behavior as introduced by [TU] . Namely, indi- 
viduals are active during seasons of fixed length T separated by winter periods. To give a 
representation of what such a system could encompass, the resource population could repre- 
sent plants producing seeds all season long, and the consumer population insects having to 
trade-off between feeding and laying eggs. All individuals die at the end of the season and 
the size of the next generation is determined by the number of offspring (seeds or eggs) pro- 
duced during the past season. The rate at which the consumers produce immature offspring 
(eggs) depends on their internal energy which can be increased by feeding. The reproduction 
of the resource simply occurs at a constant rate. 

In nature several patterns of life-history of the consumers can be singled out, but they 
almost always contains two main phases: growth phase and reproduction phase. Depending 
on initial conditions the transition between them could be strict when the consumers only 
feed at the beginning of their life and only reproduce at the end, or there could exist 
an intermediate phase between them where growth and reproduction occur simultaneously. 
Such types of behaviors are called determinate and indeterminate growth pattern respectively 

nu. 
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Time-sharing between laying eggs and feeding for the consumers is described by the vari- 
able u: u — 1 means feeding, u — on the other side means reproducing. The intermediate 
control u G (0, 1) describes a situation where, for some part of the time, the individual is 
feeding and, for the other part of the time, it is reproducing. 

At the beginning of the paper, we consider a population of consumers maximizing their 
common fitness, all consumers being individuals having the same goal function and acting 
for the common good. We will call them the residents in the following. We suppose that a 
small fraction of the consumer population may appear at the beginning of one season and 
start to behave as mutants in the main population. We study how such invasion will occur. 

If there is a large number of residents, it makes sense to assume that residents fix their 
strategy a priori and do not change it during the season. The mutant achieves better result 
than the resident using this fact and react in feedback form. Such a problem can be related 
to a hierarchical game of two players. For simplicity we study this problem in the case of a 
vanishingly small number of mutants. Such a situation corresponds to the definition of an 
evolutionarily stable strategy given by [7] , when only a small fraction of mutants is taken 
into consideration. We also investigate the fate of such a mutation in the multi-seasonal 
framework proposed by |10| . In particular we show that mutants not only can invade the 
resident consumers' population, but will also replace it in the system. Finally, we make 
some conclusions regarding the results presented in the paper. 

2 Main model 
2.1 Previous work 

At the beginning consider a system of two populations: consumers and resource without 
any mutant. As it has been stated, all seasons have fixed length T which does not change 
from one year to another (see Fig.QJ. The consumer population is determined by two state 
variables: the average energy of one individual p and the number of consumers c present 
in the system. For the description of the resource population a variable n is introduced. 
It defines the size of the population. We suppose that both populations consist of two 
parts: mature (insects/plants) and immature part (eggs/seeds). During the season, mature 
individuals can invest in immatures by laying eggs. Between seasons (at winter periods) all 
matures die and immatures become matures for the next season. 

We suppose that all consumers have arbitrarily small energy p at the beginning of the 
season. The efficiency of the reproduction is assumed to be proportional to the value of 
p: it is thus intuitively understandable that consumers should feed on the resource at the 
beginning and reproduce at the end once they have gathered enough energy. The consumer 
has a trade-off between feeding (u = 1) and laying eggs (it = 0). The variable u plays the 
role of the control. 

The within season dynamics are thus defined as follows 

p= — ap + bnu, n = —cnu (1) 
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Figure 1: Seasonal character of the behavior: populations evolve continuously during the 
season and have discrete rules for the transition from one season to another. Picture has 
been taken from [12] 

where we supposed that both populations do not suffer from mortality; a, b and c are some 
constants. After rescaling of time and state variables, the constants a and b can be eliminated 
and equations (TTJ) can be rewritten in a simplified form 



where c is represent the number of predators present in the system. 

The amount of offspring produced by individual during the season depends on the current 
size of the populations 



where consumers are maximizing the value J, the common fitness, 9 and 7 are some con- 
stants. We see that this is an optimal control problem which can be solved using the dynamic 
programming |13| or Pontryagin maximum principle |14) . Moreover, the constants c, 9 and 
7 can be omitted to compute the solution of this problem without loss of generality. 

One can also show that all the data of the formulated problem are homogeneous of 
degree one in state variables, which can be only positive numbers. This is a particular case 
of Noether's theorem in the calculus of variations about the problems whose data is invariant 
under a group of transformations [15] . Hence the dimension of phase space of the optimal 
control problem (J2K2]) can be lowered by one unit by the introduction of a new variable 
x = p/n. In this case its dynamics can be written in a form 



p = —p + nu, n 



cnu 



(2) 





x = —x(l — cu) + u 
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Figure 2: Optimal collective behavior of the residents 

and the Bellman function - a solution of an optimal control problem U(p,n,t) = J^_ t (l — 

u(s))p(s)ds with the starting point at (p(t),n(t)) = (p,ri), can be present as U(p,n,t) = 
nU(x,t). 

The solution of the optimal control problem has been obtained before [TU] and the 
optimal behavioral pattern for z = 1.5 and T = 2 is shown on Fig. [5] The region with u = 1 
is separated from the region with u = by a switching curve S and a singular arc S a such 
that 

S: x = l-e- {T ^ (4) 

S a : T-t = -lnx + — -- (5) 

xc c 

They are shown on Fig. [5] by thick red and green curves correspondingly. Along the singular 
arc S a the consumer uses intermediate control u = u: 

u = (6) 
2 + xc w 

One might identify a bang-bang control pattern for short seasons T < T\ and a bang- 
singular-bang pattern for long seasons T > T\. The value Ti is equal to 

ln(c+l) + (c-2)ln2 
Tx — (7) 

and it depends on the number of consumers present in the system. 



INRIA 



Join forces or cheat: evolutionary analysis of a consumer-resource system 



7 



The optimal value of the amount of offspring produced by individual can be computed 
using this solution. In the following we concentrate on the process of mutant invasion into 
population of consumers which uses the prescribed type of behavior given on Fig. [2] 

2.2 Consumer-mutant-resource system 

Suppose that there is a subpopulation of consumers that acts as a mutant population. They 
maximize their own part of the fitness taking into account that the main population relates 
them as kin individuals. 

Denote the fraction of the mutants with respect to the whole population of consumers by 
e and variables describing a state of the mutant and resident populations by symbols with 
subindices "to" and "r" correspondingly. Then the number of mutants and residents will be 
c m = sc and c r = (1 — e)c and the dynamics of the system can be written in a form 

p r = -p r + nu r , p m = —Pm + nu m , n — —nc [(1 - e)u r + eu m ] (8) 

similarly to ([23) - The variable u m 6 [0, 1] defines a life-time decision pattern of the mutants. 
The control u r £ [0, 1] is fixed and defined by the solution of the optimal control problem 
(121151) . 

The number of offspring for the next season is defined similarly to @: 
J r = / 9(1 - u r (t))c r p r (t) dt, J m = I 9(1 - u m (t))c m p m (t) dt, J n = / -fn(t)dt (9) 



where the mutant chooses its control u rn striving to maximize its criterion J m . 

We can see that the problem under consideration is described in terms of a two-step 
optimal control problem (or a hierarchical differential game) : on the first step we define the 
optimal behavior of the residents, on the second step we identify the optimal response of 
the mutants to this strategy. 

3 Optimal free-riding 

Since 9 and 7 are some constants, they can be omitted from the solution of the optimization 
problem J m — ¥ max. In this case the functional J m /(0c m ) can be taken instead of the 

functional J m . 

Let one introduce the Bellman function U m for the mutant population. It provides a 
solution of the Hamilton- Jacobi-Bellman (HJB) equation 



dU m 

— h max 

Ot "m 



dU m dU m 

— (-p r + nu r ) + — (-Pm + nu m ) 

OPr Op m 



cm 



•nc((l - e)u r + eu m ) + p rn (l - Um) 



= (10) 
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Introducing new variables x r = p r /n and x m = p m /n and using a transformation of the 
Bellman function in the form U m {jp r ,p m , n, t) = nU m (x r , x m , t), we can reduce the dimension 
of the problem by one using Noether's theorem. The modified HJB-equation ([TU]) takes the 
following form 

+ maxi^ [-x r (l - c((l - e)u r + eu m )) + u r ] + 

Ot u m L OX r 

dU m 

j. — \-x m (l - c((l - e)u r + £U m )) + u m ] - 

U m c((l - e)u r + eu m ) + x m (l - u m ) j = (11) 

Since the boundary conditions are defined at the terminal time it is convenient to construct 
the solution in backward time r = T — t. If we denote the components of the Bellman 
function as dU m /dx r = X r , dU m /dx m = X m and dU m /dT = v, equation (jTTJ) can be written 
as follows 

TA. = —v + max <^ A,- [— x r (l — c((l — e)u r + eu m j) + u r ] + 

A m [-x m (l - c((l - e)u r + eu m )) + u m ] - 

U m c((l - e)u r + eu m ) + x m {l - u m )\ = (12) 

where the optimal control is defined as 

u m = Heav(^4 m ), A m = d%/du m = \x r ez + A m (l + x m ec) - U m sc - x m . 

One of the efficient ways to solve the HJB-equation is to use the method of characteristics 
(see e.g. [H]). The system of characteristics for equation (JT2J) reads 

x' r = —dH/dXr = x r {l — c((l — e)u r + eu m )) — u r , 
x 'm = -9H/dX m = x m (l - c((l - e)u r + eu m )) - u m , 
X' r = dH/dx r + X r dU/dU m = -X r , 

X' m = dU/dx m + XmdU/dUm = -X m + l-u m , (13) 
v' = vdT-L/dU m = -vc((l - e)u + eu m ), 
U' m = -XrdH/dX r - X m dH/dX m + u = 

- U m c((l - e)u r + eu m ) + x m {l - u m ) 

where the prime denotes differentiation with respect to backward time: g' = dg/dr = —g. 
The terminal condition U m (x r ,x m ,T) = gives that A r (T) = X m (T) = 0. Then A m {T) < 
and u m (T) = as it could have been predicted before. 
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Figure 3: First steps of the construction of the solution 
3.1 First steps 

If we emit the characteristic field from the terminal surface t = T with u r = u m = then 

A r (T) = A m (T) = 0, [7 ro (T) = 0. 
We get the following equations for state and conjugate variables and for the Bellman function 

X r X r (X')e , %m — *^m('^)^ ? A r 0, X m 1 e , Um — ^m(l 6 ) . 

From this solution we can see that there could exist a switching surface S m : 

S m : i m = l-e-( r - t ) (14) 

such that _4 rn = on it and where the mutant is changing its control. Equation (| 14[) is 
similar to But we should take into account the fact that there is also a hypersurface 
S r , where the resident changes its control from u r = to u r = 1 independently on the 
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decision of the mutant. Hence it is important to define on which surface S r or S m the 
characteristic comes first, see Fig. EI Suppose that this is the surface S r . Since the control 
u r changes its value on S r , the HJB-equation (|12p is also changing and, as a consequence, 
the conjugate variables v 1 X r and A m could possibly have a jump in their values. Let one 
denote the incoming characteristic field (in backward time) by "— " and the outcoming field 
by "+". Consider a point of intersection of the characteristic and the surface S r with the 
coordinates (x ri ,x mi ,n). Then x ri = 1 — e~ ri and the normal vector i? to the switching 
surface is written in the form 

= VSV = (dS r /dx r ,dS r /dx m ,dS r /dx m ) T = (-1,0, 1 - x ri ) T . 

From the incoming field we have the following information about the co-state 

A, r 0, A m X ri , V X mi G — 3>mi (1 ^r*i ) • 

Since the Bellman function is continuous on the surface S r which means 

— U m U m x nil (1 e ) x mi x ri . 

The gradient VU m . has a jump in the direction of the normal vector *&: V£/+ = VJ7~ + kd. 
Here k is an unknown scalar. Then 

A + = -k, A+=£ ri , v + = x mi (l - x Tl ) + k{l — x ri ) (15) 

If we suppose that the control of the mutant will be the same u+ = (in this case 
should be negative), the HJB-equation (fT2|) has the form 

+ X+[-x ri (1 - (1 - e)c) + 1] - X+x mi (1 - (1 - e)c) - (1 - e)zU m + x mi - (16) 

By the substitution of the values from (fT5"|) to the equation (|16p we get 

k[-2{\ - x ri ) - x ri {\ - e)c] = 0, 

which leads to the fact that k = and, actually, there is no jump in conjugate variables. 
They keep the same values as (|T5|) and = 

But let one suppose that the mutant reacts on the decision of the resident and also 
changes its control on S r from u~ = to = 1. This is fulfilled if the inequality .A+ > 
holds. 

The HJB-equation (fT^]) has the form 

+ A+[-a; ri (1 - c) + 1] + A+ [-.x mi (1 - c) + 1] - zU m = . 
Substitution of the values u + , A+ and A+ from (fT5|) gives 

= 7 r 
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Figure 4: Some family of the optimal trajectories emanated from the terminal surface 



and 



{x ri ^mi ) 



(1 - e)x ri c+ (1 - x Tl ) 



x ri c + (1 



) 



which is positive when x ri > x mi . On Fig. [3] this corresponds to the points of the surface 
S r which are below the magenta line: x r — x m = 1 — e~ r . For the optimal trajectories 
which go through such points u r (ri + 0) = u m (Ti + 0) = 1. One can show that there will be 
no more switches of the control. But if we consider a trajectory going from a point above 
the magenta line then u r (ri + 0) = 1 and u m (ri + 0) = and there will be a switch of the 
control u m from zero to one (in backward time). After that there will be no more switches. 

Now consider a trajectory emitted from the terminal surface which comes to the surface 
S m rather than to the surface 5V at first. In this case the following situation as it is shown of 
Fig. 4 takes place. One might expect to have a singular arc there. Necessary conditions 
for its existence are the following 



n = o = n a + A m u r 



Hq = — v — \x r — A, ; 



(17) 
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A m — = X r x r ec + \ m (x m ec + 1) - ecU m - x m (18) 

A' m = {Am H } = Q = A m i (19) 
where the curled brackets denote the Poisson (Jacobi) brackets. If £ is a vector of state 
variables and ip is a vector of conjugate ones (in our case £ = (x r , x m , t) and ?/> = (A r , A m , f)), 
then the Poisson brackets of two functions F = F(£, %p, U m ) and G = G(£, -0, C^m) are given 
by the formula 

{F G} = (Ft + V F Um , G ) - (i^, G { + V Gc/ m ) . 
Here (•, •) denotes the scalar product and e.g. 

F^ = dF/di/j = (dF/d\ r ,dF/d\ m ,dF/dv) T . 

After some calculations the expression (fTO|) takes the form 

Ami = vec + x m + X r x r ec — (x m + 1)(1 — A r ) = (20) 

We can derive the variable v from equation (|17p and substitute it to the last equation 
P0|). We get 

^ml *^?7i 1 ~h A771 . 

This leads to A m = 1 — x m and 

+ eU m + (1 - X m )(x m £C+ 1) 



A, 



which can be obtained from equation (|18|) . 

To derive the singular control u m = u m € (0,1) along the singular arc one should write 
the second derivative 

A m = = {{A m n}H} = {{A m n }{Ho + AraUm)} = { {A»^0 }H } + {{A m H }Am}u m ■ 

Then 

{{AmHo}Ho} 2x m . . 

Um " {Am{AmU }} " 2 + X m£C 1 ' 

which has the same form as ©• 

The equation for the singular arc can be obtained from dynamic equations (|13p by 
substitution u r = and u m = u m from (|2ip : 



2 

X Yy-i £ C 



, a: m (T = ln2) = l/2. 



"771 Q | 

Finally, we have the analogous expression to ([5]) 

Sf: T - t = - in x m + — — (22) 

for e 7^ 0. If £ = the surface S m is a hyperplane x m = 1/2. 

After these steps we have the structure of the solution shown on Fig. [5] 
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Figure 5: Construction of the part of a singular arc S" 



3.2 Optimal motion along the surface S° 

Let us consider the surface S% which is shown on Fig. [5] by green color and separates the 
domain u r = 1 from u r = 0. This leads to a chattering regime with the motion along this 
surface. The solution of the dynamic equations can be understood in Fillipov sense. 

Suppose that the hypersurface S° is also divided into two regions: the region where the 
mutant uses u m = and the one where u m = 1. They are separated by a switching curve 
S and by a singular arc S" which completely belong to the surface S° . Along this surface 
the resident uses an intermediate control u r = u r resulting from the chattering regime with 
simultaneous switches from one bang control u r — 1 to another u r — and vice versa. 
We suppose that the trajectory can be forced to stay on the surface S° by the resident 
independently on the action of the mutant. This means that if we derive the control u r from 
the dynamic equation 

2 

X C 

x' r = - - y = X r (l - C((l - e)ll r + eu rn )) - U r 
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as 

2x r (l ~t- x r c) x r scu m /oQ^ 

Ur ~ (1 + (1 - e)x r c)(2 + x r c) l + (l-e)x r c 

Then u r e (0, 1) for all points belonging Sy and for all possible values u m £ [0, 1]. 

To identify for which parameters of the model this is possible, we may notice that u r as 
a function of u m is linear and decreasing. Moreover 

i — c\\ 2x r (l -\- x r c) 1 -\- x r c 

U r (U, n = Uj = (1 + (1 _ £)xrc)(2+Xrc) S ZX rJ^- c S • 

since x r < 1/2. Therefore only the condition u r > could be violated for some values of u m . 
To define the limiting value u m for which u r (u m < u m < 1) < one can use the following 
condition: u r (u m = u m ) = 0. This gives 

2 1 + x r c 



ec 2 + x r c 



If such value u m is outside of the interval (0, 1) then the condition u r £ (0, 1) holds for any 
x r belonging to S° ■ This occurs if 

e<l/c (24) 

In this paper we consider only the values of e satisfying (|2~4"|) . This has a biological 
explanation since for sufficiently large e the resident should react to the behavior of the 
mutants who does significantly significantly affect the dynamics of the system. A fixed a 
priory strategy of the resident does not make sense in that case. 

If we consider a field belonging to the surface S°, the gradient of the restriction of U m 
to that manifold is defined only in the co-tangent bundle. A safe representative requires a 
term kn s be added to the adjoint equations of the characteristic system, where n s is the 
normal to S° . The constant k should be chosen to keep the adjoint tangent to it. But we 
can notice that the surface S° does not depend on x m -coordinate. Since only the dynamics 
A m plays an important role for us and the corresponding term is equal zero, this notion can 
be neglected. 

The control u r = u r (x r , x m , r, u m ) is defined in feedback form, e.g. depends on time and 
a state of the system. The corresponding Hamiltonian (|12p is changed to 

T~L *H. (x r , x m , Ujn , A, A m , v : ii r (x r , x m , t, u m ) , u m ) . (^^) 

The coefficient multiplying the control u m is also changed to 

i dH _ A m (l + x r {\ - e)c + x m ec) - ecU m _ 
OU m 1 + (1 — ejx r c 

In this case the switching surface S can be defined by the condition A m = 0. The singular 
arc S a - by the following conditions 

H = 0, Xn = 0, A' m = {A m H} = . 
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The intermediate control u m can be obtained from the second derivative 

A; n = {{A m n}n} = o . 

We can write analytical expressions for S and S a but they look quite complicated. To make 
things simpler, let us consider first a particular case of vanishingly small values of e and 
study the optimal behavioral pattern. 

3.3 Particular case of a vanishingly small population of mutants 

We have e = and the chattering regime of the resident along the surface S° results in 
u r = u r coinciding with ([5]): 

2x r 

Ur = -Z 

2± I XyC 

which does not depend on the action of the mutant. In addition equations (|25p and (|26p 
take the following form 

" -X ^ x j,c , f 2 x ipC \ 2 ic <j- c , . , . 

ft = -V + 7— V Am ~ x m— h% - ^mTT" h X m (l - U m ) , (27) 

z + x r c \ z + x r c / z + x r c 

If the trajectory goes from the point x^ = x m (ln2) > 1/2 then u m = and the system 
of characteristics for the Hamiltonian (|25[) is written in the form 

, x r c , 2 x r c 

X-m ~ ; , A™ - A 77-1 ~r 1 . 



2 7 ~"m O i 

I XyC Z I XyC 

^iXyC TTl TT *2iXyC 



m j 



Z I XyC Z I Xf^C 

with boundary conditions 

t = 1ii2, av = 1/2, s m = < l7 A TO = l/2, v = x° m l% U m = x 2 m /2. 
Then A m = 1 — e~ T and the switching curve S has the form 

x m — 1 e , 

in addition to r = — mx r + 2/{x r c) — A/c. Thus S = S m n S° . 

The switching curve S ends at the point with coordinates (x r2 , x m . 2 , ra) where the char- 
acteristics become tangent to it and the singular arc S a appears. Before the determination 
of the coordinates of this point let one define the singular arc S a . From equations ([2"T)l and 
1(28|) we get 

j\yXy,C 2 XyC ^iXyC 

V - ~~7 j A m X m — ■ U ! "t - %m ■} A rr} - X m 

Z I XyC Z I XyC Z | XyC 
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Figure 6: Optimal behavior inside the surface S% 

along the singular arc. Substitution of p.3[) into equation A' m = gives 

2 + x r c 
x m = - . 

Alongside, the intermediate control u m can be derived from A'^ = and it is equal to 

1 



2 + x r c 

which is positive and belongs to the segment between zero and one. 

We see that the coordinates x r2 , x m2 and can be defined through the following equa- 
tions 

2 + x r2 c _ T2 2 4 
x m 2 = r = 1 - e - , r 2 = - lnav 2 H , 
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Figure 7: Optimal behavioral pattern for c = 3 

which comes from the fact that the point (x r2 , x m21 T2) belongs to S a and it is located on 
the intersection of the curves S a and S. The result is illustrated on Fig. [5] 

We can show also that the surface S° can be extended further with comparison to the 
situation on Fig. O Indeed, the following conditions are fulfilled for the region with u r = 0: 



H 

Therefore 



V X r X r X m X m -\- X m 0, >Am 'Vn X m 0, >A« m . 

=0 



and from the condition A' m = 0: — 1 + 2x m = 0, which gives x m — 1/2. 

Consider now the region with x r smaller that the ones on the green surface Sy (see 
Fig- EJ) ■ There is a switching surface which extends the surface S m and it is defined by the 
same equation (fH)) . But there could exist a singular arc starting from some points of 
S m . To check this we have to write the following conditions 



U 



= -v - \ r (x r (\ - c) - 1) - \ m x m (l - c) - zU m + x m — , (29) 

u r — l 
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*A'm A m X m 0, ^ 5 (^0) 

which give a possible candidate for a singular arc 

s°. x -_L_ 

2 — c 

We see that its appearance is possible only for c < 2. In addition, the motion along this 
surface occurs with control u m = (1 — c)/(2 — c) which also gives the restriction on the 
parameter c that c < 1. For c > 1 the structure of the solution in the domain below the 
surface S° is simpler and consists only of the switching surface S m , see Fig. [71 

3.4 Computation of the value functions in case of e = 

Without loss of generality we can assume that at the beginning of each season the average 
energy of the population of consumers is zero: x r (0) = x m (0) = 0. Therefore we should 
take into account only the trajectories coming from these zero initial conditions. The phase 
space is reduced in this case to the one shown on Fig. [5] One can see that there are three 
different regions depending on the length of the season T. If it is short enough T < T\ 
(where the value T\ has been defined in ((7J), then the behavior of the mutant coincides with 
the behavior of the resident and the main population can not be invaded: the amount of 
offspring produced by the mutant is the same as produced by the resident. If the length of 
the season is larger than T > Ti there is a period of the life-time of the resident when it 
applies the intermediate strategy and spares some amount of the resource for its future use. 
The mutant is able to use this fact and there exists a strategy of the mutant that guarantees 
better result for it. 

Let us introduce the analogue of the value function U m for the resident and denote it as 
U r : ^ 

U r {Pr,p m ,n,t) = / p r (s)(l - u r (s)) dt . 

JT-t 

The value U(0, 0, n(0), T) represents the amount of eggs laid by the resident during the season 
of length T. Its value depends on the state of the system and the following transformation 
can be done 

U r {p r ,p m ,n,t) = nll r {x r ,x m ,t) . 

In the following we omit some parameters and write the value function in the simplified 
form U r (T) = /7 r (0,0,T) where the initial conditions x r (0) = x m (0) = have been taken 
into account. 

In the region A (see Fig. [5]) the value functions for both populations of mutants and 
residents are equal to each other 

U m (T) = U r (T) = Xl e^ T - T ^ . 
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Figure 8: The reduced optimal pattern for trajectories coming from zero initial conditions 
x r (0) = x m (0) = and c = 3 
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Figure 9: Optimal free-riding of the mutant 
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Here the value t\ can be defined from the intersection of the trajectory and the switching 
curve S r n S m : 

e (c-l)(T-ri) _ i 



c-1 

To obtain the value functions in the regions B and C let one solve the system of charac- 
teristics (|3.3p in case when the characteristics move along the surface S° and u m = 1. This 
leads to the following characteristic equations for the Hamiltonian (12~T1) : 

, x r c , 2 x r c i x r c 



~r i m 9 _i_ m 9 _i_ 

We can rewrite them in the form 

dx 2 I XpC x y-fi 2 ■ x y c ^ U^ui 2 u\ 
dXf Xy,c dXf x>p 

and consequently 

£ m = Cia^e 1 " + x r z + 1, U m = Cix 2 r , Ci,C 2 = const. (31) 

where Ci and C2 are defined from the boundary conditions while the equation ^ is also 
fulfilled. 

Along the singular arc the mutant is using the intermediate strategy (|21l) . In this 
case 

+ x r c 

U m U m ZU r T" .x m ^l — U n i) — (7.^ 



2 I c 



Since x' r = — 2 ^ c - we have 



dU m 2U m (1 + x r c)(2 + x r c) 
dXf x<p 4^Z/j„c 

Then 

^2 4 + 3x r c(3 + 2x r c) ^ . 
C/ m = C 3 x 2 r + '-, C 3 = const . (32) 

We undertake now to compute the limiting season length T2 that separates the region B 
from the region C. The coordinates of the point B have been obtained before and satisfy the 
equations (|3.3p . To define the coordinates of the point (1^,1^,72") of intersection of the 
optimal trajectory with the curve AD let use the dynamics of the motion along the surface 
S° with u r = u r and u m = 1 (j3Tj) : 



x m = C\x 2 .c T + x r z + 1, C\ = const , 
where C\ should be chosen such that 

Xrri2 — C\X r e " ~\~ X r2 C ~\- 1, x m2 - 1 e 
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Then 

= (x r2 C-2)(3x r2 C + 2) 

After that the coordinates x% 2 , x a m2 and t% can be defined from the following conditions 

„. 2 .a - . - . - 2 4 



<n 2 = A = Ci « 2 ) 2 c T ^ + < 2 c + 1 , t° = - In z° 

The boundary value T2 can be obtained as 

rj , hi« 2 (c-l) + l) 



To 



c- 1 



Now compute the value functions U r (T) and U m (T) for the region B (Ti < T < T2), 
where only the mutant keeps the bang-bang type of the control. For the resident population 
we have 

U r (T) = U r2 e-< T ~ T2 \ U r2 =x r2 (l-x r2 )+ 1 ~ 2Xr2 . (33) 

c 

where the point with coordinates (x r2 ,x r2 ,T2) defines the intersection of the trajectory and 
surface S": 

2 4 e (c-i)(T~r 2 ) _ l 

T 2 = - lnx r . 2 + - , X r2 = . (34) 

C C C _L 

For the mutant population the value function U m in the region with u = u and u m = 1 
satisfies the equation coming from (pJTj) : 

te i} =<(^) 2 , (35) 

where (x ri , x mi , t\) is a point of the intersection of the trajectory with the curve AB (see 
Fig. [8]). Using (|35[) and notations of ((34)) . we can write 



U (T) = U e" c(T_T2) U - x 2 ('- 



l 2 



which is analogous to ([3"3"1) . 

For the region C the value function for the resident has the same form as in (|3"3")) but it 
has a different form for the mutant. Suppose that the optimal trajectory is coming to the 
surface S" at point with coordinates (x r2 ,x m2 ,f2). Then the Bellman function at this point 
is equal 



2 4 + 3i r2 c\ 3x r2 c(2x r2 c + 3) + 4 



c 



,16 24xf 2 c J 24x r2 c 

which is written using Q32p with definition of the constant C3 from the given boundary 
conditions. 
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Figure 10: Difference in value functions for the resident and the mutant (c = 3) 



When the optimal trajectory moving along the surface S a intersects the curve AD at 
some point with coordinates (x° , , f%) (see Fig. [5]) the Bellman function can be expressed 
as follows 

Then 

U m (T) = U^ 2 e-< T -^ . 

The difference in the values functions (amount of offspring per mature individual) of 
the mutant and optimally behaving resident is shown on Fig. 1101 In way we can derive 
the expressions for the number of offspring produced by the resource population during the 
season. 

3.5 Generalization on the small enough but non-zero values of e 

In this section we consider a case of non-zero e but such that the condition (|24[) remains 
fulfilled. This means that the trajectory coming to the singular surface S° does not cross it 
but moves along it due to the chattering regime applied by the resident (f2~3"|) . 

In this case the phase space can be also divided in two regions: the points with x r smaller 
or larger than the ones on S° . In each of this region the structure of the solution has similar 
properties as in the case considered above when e is arbitrary small. Inside the surface S° 
the optimal behavior has also a similar to a previous case structure. 

In the region with values x r larger than the ones on the surface S° there is a part of the 
switching surface S m and a singular arc Sf where the mutant uses an intermediate strategy. 
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The surface Sf can be denned through the expression (|2"2")> . In the other region we also have 
a part of S m and a singular arc S'f which is different from «Sf and could not exist for some 
values of the parameters of the problem c and e. 

To identify the values for which the surface S*^ is a part of the solution let us write 
necessary conditions similarly to (|29H30|) : 



H 



0, An = 0, A' m = {A m H} = . 



Using these equations we are able to obtain the values of A r , A m and v on the surface and 
substitute them into the second derivative A" n = {{Am^H} = to derive the expression 
for the singular control applying by the mutant on this surface: 

2x m - (1 - e)c(l + x m ) 

Urn = —Z 7Z r~ ■ (36 

2 - (1 - e)c + x m ec 

There are several conditions which should be necessarily satisfied. First of all, the control 
should be between zero and one 

< 2x '»~( 1 ~ £ ) c ( 1 + :E "') < : (- 37 ) 
— 2 — (1 — e)c + x m sc ~ 

Second of all, the Kelley condition should be also fulfilled p^l p. 200]: 

d d 2 dU 



du m dt 2 du Ti 

This leads to the inequality 



= {Am{A m n}} < . 



2-{l-s) + x m ec>0 (38) 

In particular, both conditions ([37)) and ((38)) give that x m < 2/(2 — c). 
To construct the singular arc S% we should substitute the singular control u m from (|37|) 
and u r — 1 into dynamics (1131) : 

x'm = x m(l - c((l - e)u r + eu m )) - u rn 

with boundary conditions obtained from the tangency condition of the optimal trajectory 
coming from the domain u m = u r = 1 on the switching surface S m : 

-m(i- 1 ^ 



2-c(l-£)V 2-c(l-e) 
Such tangency takes place only if 
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Figure 11: Structure of the optimal behavioral pattern for c = 1.25 and e = 0.35 

that comes from the condition that a singular surface S m exists only for < x m < 1. This 
gives the following inequality 

1 -c(l -e) > 

for the existence of the surface S% . One can check that the inequalities (|37H38|) are fulfilled 
as well. The result of the construction of the optimal pattern for some particular case is 
shown on Fig. [TTJ 

4 Long-term evolution of the system 

|10| introduced model (J2J) as the intra-seasonal part of a more complex multi-seasonal popula- 
tion dynamics model in which consumer and resources survive during one season only. They 
considered that the (immature) offspring produced by the consumers and ressources through 
some season i and defined by equations ^ would mature during the inter-season to form 
the initial consumer and resource populations of season (i+1). Up to some proportionnality 
constants accounting for the efficiency of the reproduction processes as well as overwintering 
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mortality, jTU] obtained the following relation between the number of consumers of season 
(i + 1) and the initial number of resources of season (i + 1): 



with Ji and J n s defined in equations (J3|). 

In the present mutant invasion model, the total consumer population is structured into 
(1 — £i)ci residents and EiCi mutants that have different reproduction strategies. Taking into 
account this structure and assuming that mutants' progeny is also composed of mutants, we 
get the following inter-seasonal part for the mutant invasion model. 

c n+i = oiU r (ci,£i,ni,T), c mi+1 = aU m (ci,ei,rii,T), n i+ i = f3V(ci,ei,ni,T) 

where the values U r , U m and V denote here the number of eggs/seeds produced by each 
(sub-)population: 



U r , U m and V can be computed from the solution of the optimal control problem © 
with the dynamics (0). Their values depend on the strategy chosen by the mutant and the 
resident, the length of the season T, the values c, and £j and initial conditions which are 
p r (0) = p m (0) = and n(0) = n,. For a particular case e = the values U r and U m were 
derived analytically in subsection !3.4l In the following we investigate numerically this model 
on typical example; in particular we are interested in the long term fate of the resident and 
mutant consumer populations. 

A previous investigation of the inter-seasonal model with collective optimal behavior of 
the consumers (i.e. there are no mutants such that e = 0) has shown that the behavior 
of the system in long-term perspective could have very rich properties. Depending on the 
parameters of the problem, the value of j3 and the length of season T, there could be an 
extinction of the resource or a or blowing up of the system (which leads to the suicide 
of the consumers). The system could also tend to some stable periodic behavior or to a 
globally asymptotic equilibrium. The last two cases illustrate a possible co-existence of the 
interacting species |10j . 

Here, we follow an adaptive dynamics like approach and consider that the resident con- 
sumer and the ressource population are at a (globally stable) equilibrium and investigate 
what happens when a small fraction of mutants appear in the resident consumer population. 
We actually assume that resident consumers are "naive" in the sense that even if the mutant 
population becomes large through the season-to-season reproduction process, the resident 
consumers keep their collective optimal strategy and take mutants as cooperators, even if 
they do not cooperate. 



Cj+i = oiJi, n i+ i(t = 0) = (3J n ,i, 
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Figure 12: Fraction of the mutants in the population of consumers 
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Figure 13: Effect of the mutant invasion on the system 
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We investigated numerically a case when a — 2, (3 = 0.5 and T = 4. The system is 
near the long-term stable equilibrium point c = 0.9055 and n = 1.0848 as at the beginning 
of some season a mutant population of small size c rn = 0.001 appears. The mutant pop- 
ulation increases its frequency within the consumer population (see Fig. I12p and modifies 
the dynamics of the system (Fig. HU]) . The naive behavior of the consumers is detrimental 
to their progeny: along the seasons, mutant consumers progressively take the place of the 
collectively optimal residents and even replace them in the long run (Fig. 112(1 . making the 
mutation successful. We should however point out that the mutants' strategy as described 
in ((9]) is also a kind of "collective" optimum: in some sense, it is assumed that the mutants 
cooperate with the other mutants. If the course of evolution drives the resident population 
to and only mutants survive in the long run, this means that the former mutants become 
the new residents, with actually the exact same strategy as the one of the former residents 
they took the place from. Hence they are also prone to being invaded by non-cooperating 
mutants. The evolutionary dynamics of this naive resident-selfish mutant-resource appears 
thus to be a never-ending process: selfish mutants can invade and replace collective optimal 
consumers, but in the end transforms into collective optimal consumers as well, and a new 
selfish mutant invasion can start again. We are actually not in a "Red Queen Dynamics" 
context since we focused on the evolution of one species only, and not co-evolution [T7]. Yet, 
what the Red Queen said to Alice seems to fit very well the situation we just described: 
"here, you see, it takes all the running you can do to keep in the same place" |18) . 
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